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$_i i I. INTRODUCTION 

Oh; 

The time-dependent quantum mechanics formalism has been extensively used in the recent years in several problems 
of physics and chemistry. Since the work of McCullough and Wyatt jij showing the interest of dealing directly with 
non-stationary states in calculating reactions probabilities, a consequent number of works has been devoted to the 
wave-packet dynamics. The reader could find the more recent developments in the reviews |^-^] and references therein. 

When solving the time-dependent Schrodinger equation in a finite volume, one is faced with the problem of im- 
plementing the proper boundary conditions at finite distance. These boundary conditions are usually those of a 
free wave outgoing propagation. For the stationary solutions they are given by fixing the value of the logarithmic 
normal derivative on the boundaries of the interaction domain. However for the case of non-stationary wave packets 
' a satisfactory solution is not known. 

The usual way to overcome this difficulty is to impose the nullity of the wave on the boundaries and push them 
far enough to avoid the effect of parasite reflections in the region of interest H . Another way, widely used in other 
branches of physics and mechanics as well, is by means of the so-called absorbing boundary conditions [@-f|. They 
^ c~| indeed minimize parasit reflections but perturbe the dynamics near the boundary and do not provide a solution of 
<*S i| the initial equation in its neighborhood. 

A third way is by the so called wave function splitting algorithm which consists in removing from the total wave- 
function the part localized outside the interaction region before it reaches the boundary [|l0| . Other possibilities based 
on the interaction picture formalism have also been developed . 

Although the above quoted methods provide a practical way to solve the problem they do not give a theoretical 
solution of it. The aim of this paper is to formulate an exact boundary condition (EBC) at finite distance for 
outgoing solutions of the time-dependent Schrodinger equation. Consequently, it would allow a solution of time- 
dependent quantum mechanical problems in a finite spatial domain without any approximation in the boundaries. 
This theoretical result is completed by implementing the derived boundary conditions in some practical calculations 
concerning the more usual problems of the wavepacket dynamics. 

The paper is organized as follows: in Section 2 we derive the general expression for the EBC which has to be 
imposed at the boundary of an arbitrary integration domain containing the interaction. This conditions turns to 
be non local on time and allows a solution limited to the interaction domain of the time-dependent Schrodinger 
equation. In Section 3 the boundary conditions are reformulated for a discrete-time evolution problem in the frame 
of the Crank-Nicholson scheme. Section 4 is devoted to the one dimensional Schrodinger equation for which the use 
of these conditions is particularly simple. Several illustrative examples are given in Section 5: the spreading of a free 
wave packet, the scattering of a wavepacket on a static potential, the behavior of an initially localized wave packet 
under a time-dependent perturbation. In the last example we derive explicit solutions of the Schrodinger equation 
in a time-dependent delta potential. The applicability of the formulated conditions in a two-dimensional problem is 
shown in Section 6. A discussion on the interest and limitations of the EBC is finally presented in the conclusion. 

II. BOUNDARY CONDITIONS AT FINITE DISTANCE FOR THE TIME-DEPENDENT SCHRODINGER 

EQUATION 

Let us consider the dimcnsionless time-dependent Schrodinger equation in a d-dimcnsional configuration space 
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[id t +V%-V{x,t)]V{x,t) = (1) 

The potential V , eventually time-dependent, is assumed to vanish outside a finite domain Di named the interaction 
domain. 

We wish to obtain a solution of ([!]) by solving the equation in a finite spatial region Di, called hereafter the 
integration domain. The domain Di has to be chosen such that it contains Di, i.e. such that its complementary D e 
with respect to the total configuration space D is free of interaction. We have thus D — Di[j D e and Dj C Z)j ( see 
figure Q). 

The solution of (|lj) in Di requires to formulate and implement the proper boundary conditions on its border, a 
surface E. For this purpose let us first consider the solution of (|l]) in D e . The domain D e is, by definition, free of 
interaction and \& satisfies the free Schrodinger equation 

(id t + vl)*(x-,t) = o (2) 

Let Kq be the free retarded propagator, i.e. a solution of 

{-id t + V 2 s )K+(xb -x,to-t) = i6(xo - x)S(t - t) (3) 

with x G D e and (xoi^o) arbitrarily chosen. 

We have in mind the outgoing solution at (xo,io) of a time-dependent problem with initial value given by ^>(x,ti) 
for ti < to. By multiplying equation (Q) by Kq, equation (^) by \f and substracting one gets the relation 

9(2, t)S(x - x)6(to -t) = -d t [K+(xo - x, t - t)V(x, t)] 
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i^(x,t)^ s K+(x -x,to-t) (4) 



with / V g = f(Vg) - (Vf)g. 

This relation is to be integrated over the d + 1 dimensional space-time volume D e x [£j,+oo]. We will assume 
that the domain D e is connected. This assumption is justified for most of the cases except for d = 1, which will be 
discussed later. The integral of the right hand side is then transformed, by using the Green theorem, into a flux term 
through the boundary surface S. We denote by n the normal to £ pointing towards D e (see figure Q). 

Taking into account that Kq (x~o — x, to — t) — for t> to and that there is no contribution from the flux at the 
infinity for purely outgoing waves (for they have the same asymptotic behavior as the retarded free propagator) we 
obtain the following relations depending on the relative position of xo- 

If x € D e one has 



#(x"o,to)=/ K+(x - x,t -U)^{x,ti)dx (5) 
+i J dt J ^(x , t) V K^(xo - x, t - t).ndS 

The solution in the free region is thus obtained as a sum of a first term, representing the free propagation of the wave 
which was initially in the non interacting region D e , and a second term made only with the values of the wavefunction 
and its normal derivative at the boundary S at the past times. Equation (||) can be consequently used to propagate 
outside Di a solution obtained in Di and can be viewed as a generalization of the Huyghens principle in a space-time 
surface. 

If xo G E one has 



1 ^(xo,to)= I K+{xo-x,to-ti)^{x,ti)dx (6) 

-De 
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+i dt J *(x, t) V KfKxo - x,t - t).ndS 



it, Js 

This result can be obtained either as a limiting case of the preceding one @ or by directly integrating equation (^) 
for x G E. Equation (g) is the boundary condition we were looking for. It provides an integral relation between the 
value of the wavefunction on a point of E at a given time and the values of the wavefunction and their derivatives 
in E at the preceding times. Starting from a known initial state at time U, equation (^|) gives in its differential form 
the boundary condition on E required for determining at time ti + At the solution of the time-dependent Schrodinger 
equation in the finite domain Di. 
Several comments arc in order: 
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The relation (g) is non local in both space and time. It generalizes the outgoing boundary condition for the 
monochromatic free plane waves, which in the one-dimensional case and at x = a reads: 

tf(a,to) = -rrd x ^(a,t ) (7) 

%K 



• The choice of Kq ensures that only free outgoing waves will propagate across E, in the same way as a solution 
obtained in an infinite domain. This result is exact and consequently will not generate any reflection on the 
boundary. 

• The first term of its right hand side corresponds to the free propagation of the initial state which was initially 
in the external region D e . Its numerical calculation needs some care due to the oscillatory character of the 
propagator and several ways to overcome it are proposed below: 

— It is always possible to eliminate this contribution by choosing the initial time ti and Di in such a way that 
the initial state is fully in Di. 

— In the case when the initial state is very extended, it is more suitable to include it entirely in D e . Indeed the 
volume integral term gives then the value at the boundary of a free evolution of the initial state. This can 
be easily calculated in the momentum space and is known analytically for some cases of physical interest 
like gaussian wavepackets. 

— In the case when the initial state is entirely in D e , there exists another possibility to remove the volume 
term by splitting the solution into its non perturbed and scattered parts, i.e. 

* = $ + x 

in which $ is the solution of (g) which coincides at t — ti with the initial state and \ a solution of the 
inhomogeneous equation 

[id t +Vl] X &t) = V(a!,t)$(x,t) (8) 

The unknown function x vanishes at t = U and satisfies in the free domain D e the same equations than 
It will consequently satisfy on S the boundary conditions @ without the volume term. 

The preceding relations, completed by the case xq £ Di, are summarized in the following equation 



e*(aT , t ) = / K+{xh - x,t - U)^(x,U)dx (9) 
+i [ dt { ^(x,€yv *K£(xo- x,t -t).ndS 



with e defined by 



1 if a?o £ D e 
e (^o) = { if x G A 
| if x e X 

The above derivation has been done under the assumption that the potential V vanishes outside a finite domain 
Di. This condition can be however weakened in some cases. For instance in the case of an interaction made of one 
short range part V s plus a static long range part Vj, (e.g. Coulomb, polarisation or centrifugal potentials) the same 
derivation holds. The interaction domain is then defined by V s whereas D e can contain Vl all across. The only 
difference in equation (||) consists in replacing the free propagator Kq by the corresponding propagator of the long 
range potential Vl- 

In the one dimensional case the interaction domain is Di = [—a, a] and the free exterior domains D e has two 
connected components [— oo, — a] and [a,oo]. The boundary condition (^), based on the Green theorem, has to be 
applied to each of them. Let us take for instance Df = [a, oo]. At xq = a, the normal derivative of the free propagator 
on the boundary vanishes, as it can be explicitely checked on its expression (h = 2m =1): 

K+{a-x,t Q -t) =Q(t -t) 1 e^ 1 ^- (10) 

y/4lir{t - t) 
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The EBC becomes simply 



i f°° r / i 

-*(o,tb) = y K+(a-x,to-ti)y(x,ti)dx- J dtJ _ d x V(x,t)\ a (11) 

It is worth noticing that in case f j = — oo this result may be obtained in a quite straightforward way. Indeed let us 
consider the Fourier energy components of the wavefunction 

r+oa 

V(x,to)= dw*(x,w)e- lu *° (12) 



and let us impose to each of them the plane wave outgoing boundary condition (Q) on the form 

iy/uV(a,to) = d x ^{a,t ) 

Inserting this relation in (111) one gets 

r+oo /.+00 e -i u (t -t) 

*(o,to)=/ dt du d x V(a,t) (13) 

We recognize the Fourier transform of the free Green function for zero argument 

/+oo -, 
dt e-^K+ix, t) = -^e^M (14) 

and we obtain the following relation: 



*(«./„! = - / ° dtJ 1 -r d x *(a,t) (15) 
, V Tr(t - t) 



in agreement with the second term in the right hand side of (|1 1|). The first term does not contribute in the limit 
ti — > — oo due to the fact that the free propagator tends to zero. 

The boundary condition (||) generalizes the KKR |lj] method which provides a relation between the solution and 
its normal derivative at the boundary of a free domain. This method has been derived and successfully applied in the 
case of static quantum billards. The formulation we have developed in this section appears as a natural extension of 
this method for space-time billards. 

The condition ([)]) is not directly useful in a numerical solution of the problem. Indeed the singularity of the 
propagator at t = to would generate instabilities when the dependent Schrodinger equation is discretized. To overcome 
this difficulty one has to reformulate the boundary conditions in the particular approximate scheme chosen. The aim 
of the following section is to obtain such a formulation in the frame of the Crank-Nicholson scheme. 



III. BOUNDARY CONDITIONS IN THE CRANK-NICHOLSON SCHEME 



The Cranck-Nicholson |13| method for solving a discrete-time evolution problem consists in approximating the 
infinitesimal evolution operator 

|*(t + At)) = e -i # At |¥(i)) 

by the expression: 

(l + ^A*)|* m ) = (l-^At)|* TO _ 1 ) (16) 

where |^ m ) and |^I' m _i) are the wavefunctions at two consecutive times separated by At and the hamiltonian H is 
evaluated at the mean time t + 

By introducing the imaginary parameter /i 2 = ^Kt^ ecma tion ( |l6| ) can be rewritten in the form 

(fi 2 - H)V m = (fi 2 + ff)* m _! (17) 



4 



according to what the solution at time-step m is obtained by solving an inhomogeneous stationary Schrodinger 
equation with complex energy [i 2 . The source term is provided by the solution at the preceding time. 

By means of ( |l7| ) the equivalent of equation (||) for the outgoing wavefunction in the free domain D e becomes 







and the discrete free retarded propagator K p (xq — x) is defined as the solution of 

fi 2 (K p+1 - K p ) = -V%(K V+1 + K p ) 



(18) 



(19) 



normalised to Ko(x~o — x) = S(xq — x) 

For xq (fc S one obtains the analogous to the continuum case by simply replacing the time integral by a sum over 
the discrete index p 



K n {x" - x)^ (x)dV 



(20) 



p=0 



2~2$3 /_ (K p (x - x ) + K p+1 (x - X )) V (*n-p(x) + *n_p-l(x)) 



ridS 



In calculating the limit xq — > £ and like in the continuum case there is a singularity in the normal derivative of the 
free propagator. Its contribution to the surface term in the r.h.s of (|2C| ) is 

1 n_1 1 

- £(-l)*(* n _ p ( a r ) + ^-p-x^o)) = - (* n (aS6) + (-ir* (i&)) 

p=0 

for a?o S I? e and the opposite for x*o £ Di. The volume contribution can be in its turn re- written by means of the 
regularized propagator K p defined as 



K P (X) = (-l) p 5(X) + K P (X) 



and one finally arrives to 
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K n (x - x)tp (x)dV 
1 " 1 f 

Y^Yl / (k p (x-x )+K p+ i(x-x )) V (*„_p(f) + * n -p_i(x)) 



(21) 



(22) 



fidS 



p=o • 



This expression is the equivalent of (Q) for the discrete time evolution. 

To make use of the discrete EBC condition fl2^) it remains to derive an expression for the discrete-time free 
propagator K p . In the momentum space equation (|19]) reads as 



which leads to: 



K p+1 (k) 



K p (k) = 



fi 2 + k 2 
[i 2 - k 2 



K p (k) 



2\P 



(23) 



One can see in equation (J23J) that in the limit k — > oo the Fourier transform of the propagator contains the singular 
delta function manifested in (pl|). Nevertheless this term cancels in the sum of two successive propagators. This sum 
is given by the non singular integral: 



K p (X) + K p+1 (X) 



2fi 2 



■„2\ P 



/i 2 - k 2 



iiix dk 



(2TT) d 



(24) 



This integral can in principle be calculated by the contour Cauchy method. However, its value can also be obtained 
by means of the discrete Fourier tranform. We note indeed that a N discrete-time Fourier transfom with respect to 
the index p of equation (|3) gives the sum of a geometric series. With the inverse Fourier transform we obtain: 
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K p (x)+K p+1 (x) = ±'£ {^Go(tf,x) (25) 
1=0 

where <pi = 2s + i-q and -q is the required convergence factor of the geometric series which has to be chosen such that 
e ~ N v << Q+ i s the standard free Green function with energy kf = /•t 2 y^i?j-- 

Although the formalism has been developed in view of solving the Schrodinger equation in a finite domain Di, 
equation ( p(i|) with e = 1 provides the wavefunction in all the configuration space. However one can obtain a closed 
form for some observables in terms of the internal solution only. Let us consider for example the overlapping between 
two outgoing solutions $ et VP of (|l^) in the exterior domain. Standard algebra gives for the exterior part of the 
integral the expression 

; $| *)« _ ( $| tf^-i = J_ f ($„ + V(*„ + ^ n -i)ndS (26) 



which can be calculated from the solution in Di. By adding the corresponding equations for consecutive time steps 
one obt ains the desirate scalar product. This expression is used e.g. in calculating the auto-correlation function (see 



section 



. 1VE ). The case $„ = ^> n provides the rate of variation of the probability in the exterior domain by time step. 
The contribution from Di can be calculated by numerical integration and we thus have the possibility to verify the 
conservation of the total norm . 

In summary the solution of the time dependent Schrodinger equation in the interaction domain Di can be obtained 
by solving, at each time step, a stationary complex and inhomogeneous Schrodinger equation ([l6]) in Di with the 
boundary conditions given by (p2|). The local character of this equation is preserved by these conditions. We are thus 
let with the solution of a banded linear system, whatever the discretisation algorithm used for the spatial variables. 
For static hamiltonians this linear system is in addition the same at any time-step. 

IV. EXAMPLES IN THE ONE DIMENSIONAL CASE 

In the one dimensional problems, simplifications arise. Like in the continuum case, the normal derivatives of the 
propagator at the boundaries vanish. The sum of two successive propagators for xo = a can be explicitly calculated 
and gives: 

2fi 2 U 2 + k 2 \ p dk 
p?-k 2 \n 2 -k 2 ) 2tt 



K p (0) + K p+1 (0) = I -^- 2 ) - (27) 



This integral can be performed by standard methods: it vanishes for odd values of p and for p = 2q it gives 

K 2q (0) + K 2q+1 (0) = -inC q (28) 



with 



a- <29)! 



Jq (2V) 2 



The boundary conditions ( p2| ) give at x$ = a: 

f°° i 
*„(a) = 2 K n (x - a)^ (x)dx - - ^ C q d x [*„_a,(o) + * n _ a ^i(o)] (29) 

Ja V q =o 

At xq = —a one obtains in a similar way 

/—a • t 2 ] 

K n (x + a)<$ (x)dx + C i dx [V>n-2<?(-a) + * n _ 2g -i(-a)] (30) 

For illustration we have detailed some simple examples concerning the time evolution of wavepackets. The initial 
state has been taken for simplicity of gaussian form 



G 



* (ar) = 1 e l ^ {x - Xa) e (31) 



The examples have been solved by scaling the integration domain Di to the interval [— 1, +1]. The discrete-time 
Schrodinger equation ([if]) is there solved at each of the N time-steps by finite differences method with N x equally 
spaced grid points Xi. The values of the wave function at Xi are the unknowns. The required N x equations are 
obtained by validating ( |l6| ) at each grid point. However to evaluate the second derivative at the end of the grid 
the values of the wave function at the nearest points exterior of this domain is needed. These are written in terms 
of the internal points by using EBC conditions ( p9|]3C| ) with a symmetric formula to evaluate the derivative term. 
Substituting this exterior values in (Il6|) we end with a tridiagonal linear system of dimension A^. 



A. Free propagation of a gaussian wave packet 

The first example concerns the free propagation of a gaussian wave-packet. Its main interest lies in the fact that an 
analytical solution is known which allows to check the validity of our boundary conditions. In this case the integration 
domain is not fixed by the interaction but chosen such that it fully contains the initial wavepacket. 

The wave packet ( |3l| ) is initially centered (xq — 0) and, in absence of interaction, the probability density at time t 
is known to be given by : 

p{ X ,t) =| |2 = _L_ e - 1 ^ 



with <r(t) = <7o V 1 + t 2 and i = -^jxt a dimensionless time variable. If the total evolution time in this units is T, one 
has n 2 = ~ T^f- W e have taken for the initial width the value ct = 0.2. 

In figure (||) are displayed the density probabilities each 4 time steps. The time T is equal to 4, the number of time 
steps is N = 40 and iV x = 201. The wave packet velocity v is chosen such that during the time T a classical particle 
would travel half the interval [—1, +1]. The initial state is shown by dashed line. One can see by simple inspection the 
absence of any parasite reflection nor anomalous behavior. The curve corresponding to the time t = T is compared 
to the analytical result. The difference is not visible by eyes except at the vicinity of x = 1 in which it reaches its 
maximum value, less than 1%. This difference has been arbitrarily reduced by increasing the number of time and 
space grid points showing the validity of this approach. The parameters chosen in this simple example may serve to 
illustrate the efficiency and accuracy of this method. 

In the example considered above and due to the value of v, the probability flux leaves the interaction domain mostly 
at x = 1. However, due to the spreading of the wavepacket, the wave function at x = —1 is also being populated 
although its value remains very small and is not visible in the figure. By taking v=0 one has a symmetric situation 
for which the same quality of results has been obtained. 

A last check has been done concerning the total probability conservation. The contribution coming from the interior 
domain has been obtained by integrating (trapezoidal rule) the calculated wavefunction. The contribution from the 
external part has been evaluated by equation (|2^) which becomes in that case 

r + oo i 

/ (| ^ n (x) | 2 - | tf n _i(aO | 2 ) dx = T-^i^n + *„-l)*<9 x (*n + U=l (32) 



In this example the probability is found to be conserved better than 10" 
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B. Scattering by a static potential 

In the second example, the initial wave packet (|l|) scatters on an attractive potential: 

V(x) = V e~& 

with b = 0.05 and Vb = —150. The wavepacket parameters arc xo = —0.3, ero = 0.15, v = 0.37. 

The probability densities at different times are displayed in figure (||). Its initial value is plotted in dashed line 
and the attractive well in bold. The reflected and transmitted wave packets are clearly separated. The integrated 



7 



probabilites in D e and Di have been calculated respectively by (|32|) and trapezoidal rule from the solution in Di. They 
are shown as a function of time in figure (^). The contribution from the x < — 1 and x < 1 regions are respectively 
plotted in curves (a) and (b). They correspond to the reflected and transmitted wavepackets. The difference between 
both curves corresponds to the integrated probability in the integration domain. It tends to zero very slowly. 



C. Localized state under a time-dependent perturbation 

The third example concerns the evolution of a wave packet initially localized on an attractive gaussian potential 
whose amplitude varies periodically: 

V(x,t) = V (l + sm2nLut) e"! 2 " (33) 

with b = 0.05, uj = 0.05 and V = -200. 

The time of evolution is T=80, with 800 points in space and 800 time-steps. The initial state is a superposition 
of the bound states plus a small component (i=s 3%) of continuum states both in the corresponding static potential 
at t = 0. Curves (a) and (b) on figure (|5|) represent respectively the probability density in the intervals x < — 1 and 
x < 1. The difference between these curves, i.e. the probability density in the integration domain, tends to zero 
showing that the bound state is pulled out of the potential by the time dependent perturbation. The evolution in the 
corresponding static potential, shown in dashed lines, tends to a constant: the projection of the initial wavepacket 
into bound states. 



D. Tunneling 

The following example illustrates the tunnel effect through a double well. A gaussian wavepacket with uo = 0.12, 
initially centered and at rest v = evolves in a double repulsive well of the form 



V = V 



+ e 



(g + °o) 



(34) 



with b — 0.05 ao = 0.5 and Vq — 150. We have displayed in figure (||) the probability density at different times. 
The initial state is drawn in dashed line and the potential (arbitrary units) in bold. The wavefunction spreads and 
oscillates inside the potential except for a small part that tunnels through. The time evolution of the integrated 
probability densities corresponding to figure (|^) is shown in figure (Q). Curve (a) corresponds to the domain x < — 1 
and curve (b) to i < 1. The distance between both curves represents the integrated probability in the integration 
domain Di. The small oscillations correspond to the back and forth reflections of the wavepacket inside the well. 



E. Time-dependent Delta potential 

An interesting example is provided by the limiting case of a time-dependent delta potential 

V(x) = -X(t)S(x) 

In this case the interaction domain is reduced to a point and it is possible to obtain closed analytical solutions for the 
discrete-time evolution problems. 

We will consider the evolution of a state W which initially coincides with the bound state $o of the potential with 

A 2 

strength Aq. Its energy, in current units (h — 2m — 1), is Eo = — ujq = — j- and the normalized wavefunction 




By integrating the discrete-time Schrodinger equation ([17]) on both sides of the singularity one gets the discontinuity 
of the wave function derivative at x = 

A« {*n(0) + *n-l(0)}) + [9x(*n + *n-l)]S- = (35) 
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where A„ is the mean intensity of the potential between the n and n — 1 time-steps. 

The solution of (|l7| ) ^„ will be written by splitting its non perturbed $ n and perturbed x« parts 

* n = $n + Xn (36) 

By definition, $ n is a stationary solution of equation ( |l7|) with A = Ao- In the continuum case the time evolution of 
$ n would be simply given by 

$„ = e" m " oA *$ 
In the discrete-time evolution it is easy to show that 

$„ = e-™ ;e $ 

with 

e -i0 = V 2 -^o 
/i 2 + u 

The perturbated part \ n is a solution of (|l^) on both sides of x — and satisfies the boundary conditions <\2 
The sum of these two conditions gives 



2x™(0) = -- C q [d x ( Xn -2 q + X„- 2g -i)]o- (37) 
Combining this last result with equations (p5[) and (pq) one obtains at x=0 



2Xn = - C 1 [ A «(Xn-2g + Xn-2q-l) + (K - A )($„^2 9 + $ n -2 9 -l)] (38) 

" 9=0 

This is an inhomogeneous linear system for the perturbed wave function \n with the inhomogeneous term given by 
By isolating the q=0 contribution equation ( |38| ) results into an explicit recurrence relation 

n-l n-1 
2 2 

(-2ijU - A„)x„ = A 

nXn— 1 + ^ CgA„(Xn-2g + Xn-2q-l ) + ^C,(A„-A )($ n—2q + $n-2g-l) (39) 

The perturbed wave function \n at the origin can be thus easily obtained. The values at x ^ can be calculated from 
equations ( p0| ) and ( p5| ) This example illustrates how we can extract all the physical quantities from the knowledge 
of the wave function and its derivative on the boundary, even in a limiting case when it is reduced to a point. 
In figure (||) are displayed the results corresponding to a periodically driven delta potential with 

A(t) = A + ^{l-cos(0.7^ t)} 

and where loq is the Bohr frequency of the bound state for A = Ao- Curve (a) shows the square modulus of the auto- 
correlation function < $(f)|4 r (f) > and curve (b) the square modulus of the wave function at the origin normalized to 
one at the initial time. For the long times the wave function near the origin becomes proportional to the bound state 
and these two quantites are close to each other. There are 1000 time steps for 40 perturbative pulses. We note that 
after a transient period, the wave function is near of an eigenvector of the Floquet operator jl4| with complex energy. 
The periodically driven delta potential is used as a model of weakly bound atom (ID delta atom) when studying 



the tunnel effect in a static exterior field |15|. In this case the exact propagator is known explicitly. In our case there 
is also an explicit solution. 
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V. A TWO DIMENSIONAL EXAMPLE 



We present in this section a straightforward application of the formulated boundary conditions for a two dimensional 
problem. We will illustrate it by considering the free evolution of a two dimensional gaussian wave packet crossing 
obliquely the x = 1 axis. The wave packet has initial widths uq x — <7o y — 0.2 and is initially centered at the point 
(0, 1) of the XY plane with velocity v y /v x — 3/2. We push to the infinity the domain for the y variable and consider 
as integration domain the band Di = [—1, +1] x R 

Let us apply the boundary conditions ( |22] ) to a point xq in the line 7 a = {(a,y) € R 2 : y G (— oo,+oo)}. The 
integral term has only contributions coming from 7 a itself and is thus reduced to a one- dimensional integral although 
with a two- dimensional propagator of non- vanishing argument K p (X). 

1 1 v— I f+°° 

~V n (a,y ) = Y^Yl / dy [(K p (x-x a )+K p+1 (x-xo))V{^n-p(^) + ^n-p-i(x))} (40) 

" p=0 "'~ oc 

with x = (a,y) and xq = (a,yo). 

An efficient way to solve (|40|) is by Fourier transforming the wavefunction with respect to the y variable 

/+oo 
dk v e lk « y V n (a,k y ) 
-oo 

and impose the desired boundary condition to each of its Fourier components ^jn(a, k y ). The problem is then formally 
equivalent to a one- dimensional one with a propagator given by expression ( |2^ ) where the k 2 variable is replaced by 
k 2 + k 2 . The resulting integral is no longer analytic as in the one dimensional case but can be easily evaluated by 
means of a discrete Fourier tranform according to equation ( ps|) . 

The numerical results are displayed in figure ^[ They correspond to the solution of the discrete time Schrodinger 
equation in the region x S [—1,-1-1] and y G [0,5] with N = 100 time steps. The number of points in x direction is 
N x = 101 and N y — 45 points were required for an accurate Fourier transform of \t n (a, y) amplitudes. The countour 
plots of the probability density are shown for six equidistant times. 

The figure shows the spreading of the wave packet in both directions and the crossing of the x = 1 boundary 
without any parasit reflections. The comparison of calculated values with the exact analytic expressions showed a 
relative accuracy at the maximum of the wave packet better than 1%. 

VI. CONCLUSION 

Exact boundary conditions at finite distance for the time-dependent Schrodinger equation have been derived. These 
conditions allow a resolution of time-dependent quantum mechanical problems in a finite spatial domain containing 
the interaction without introducing any parasite reflections, nor changing the dynamics near the boundary. They are 
specially useful to observe long time evolution of waves packets in a finite space domain. 

The value of the wavefunction outside the integration domain can be easily obtained from the values of the wave 
function and its normal derivative at the boundaries for all the anterior times. 

These conditions have been explicitely reformulated in a time-discrete dynamics and implemented in a numerical 
algorithm providing a numerical solution of the more usual one-dimensional problems. They result in solving at 
each discrete time step a stationary Schrodinger equation with complex energy in a finite spatial domain. The local 
character of this equation is preserved by the boundary conditions and we are thus let with a band matrix problem 
whatever the discretisation algorithm used. 

The economy resulting when solving the problem in a reduced spatial domain allows to considerably increase the 
number of grid points and to use for a same accuracy less powerful but simpler algorithms. The numerical results 
obtained with the Cranck Nicholson scheme arc found to be very satisfactory. Other time and space discretisation 
algorithms can however be used as well at the cost of the simplicity. 

The derived boundary conditions can be applied to more general situations than those considered in the examples. 
The case of higher dimensionality and/or of coupled equations, for instance, can be treated without modifying the 
formalism. A numerical application of the free evolution of a two-dimensional gaussian wave packet is given. Long 
range static potentials may be taken into account as far as an analytic solution is known for the corresponding 
propagators. The generalization to the time-dependent treatment of the three body problem is straightforward and 
could be of some interest in some molecular physics problems ]l7| or a doorway to solve the still open question of 
Coulomb problem above threshold. 



10 



The interest in formulating exact boundary conditions at finite distance for time-dependent problems recently 
appeared in the field of statistical mechanics and conclusions similar to those presented in this work have been reached 
in But this formulation is also relevant in other branches of physics than those governed by the Schrodinger 

equation. In particular for the diffusion equation, for which the same results hold with imaginary time. The case of 
hyperbolic equation of interest in many fields of physics can be treated in the same footing Kg]. 
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FIG. 1. Space-time representation of different equal-time sections of the configuration space D. The shadowed region repre- 
sents the interaction domain Di . The time-dependent Schrodinger equation is solved in the finite domain Di 
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FIG. 2. Free evolution of a gaussian wave packet. It has been obtained by solving the time-dependent Schrodinger equation in 
a limited region of the configuration space with the boundary conditions given by (^). Dashed curve represents the probability 
density of the initial wavepacket for <jo = 0.2. The time evolution is plotted at ten different times showing the spreading of the 
wavepacket and its global displacement without any parasite reflections at the boundaries x = ±1 
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FIG. 3. Scattering of a gaussian wavepacket on an attractive gaussian potential drawn in bold line. The initial state, in 
dashed line, is centered at xq = —0.3 and has ao = 0.15 and v = 0.37. It splits into a transmitted and a reflected part which 
goes accross the boundaries of the integration domain without any reflection. 
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FIG. 4. Time evolution of the integrated probability densities for the scattered wavepacket of fi gure (|3|). Curves (a) and 
(b) represent the probability density in the intervals x < ~ 1 and x < +1, calculated by equation (|32[). They tend towards 
the corresponding transmission and reflection coefficients. This convergence is relatively slow due to the slow spreading of the 
initial state in the integration domain Di — [—!,+!]. 
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FIG. 5. Evolution of a localized state in the time-dependent potential (p3|). The initial state, a gaussian wavepacket with 
<7o = 0.10, xo = 0,v = 0, is a superposition of bound states plus a small contribution from the continuum, both in the 
corresponding static potential at t = 0. Curves (a) and (b) represent respectively the probability density in the intervals 
x < —1 and x < 1. The difference between these curves, i.e. the probability density in the integration domain, tends to 
zero showing that the bound state is pulled out of the potential by the time dependent perturbation. The evolution in the 
correponding static potential is shown in dashed lines. 
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FIG. 6. Tunneling through a double repulsive well (bold curve). The wavepacket is initially centered (dashed line) with 
(To = 0.12 and v = 0. It spreads and oscillates inside the potential except for a small part that tunnels through 
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FIG. 7. Time evolution of the integrated probability densities corresponding to Fig. (g) are shown. Curve (a) is the integrated 
probability in the region x < — 1 and curve (b) in a; < 1. The distance between both curves corresponds to the integrated 
probability in the integration domain Di. The small oscillations are due to the back and forth reflections of the wavepacket 
inside the well. 
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FIG. 8. Time-evolution of a bound state in a pulsed delta potential V{x) = — A(t)5(x). Solid curve represents the probability 
density at the origin. Dashed curve represents the squared modulus of the autocorrelation function |(*o| *n)| 2 - Abscisse units 
are the number of time-steps. The total time evolution corresponds to 32 Bohr periods. The frequency of the potential pulse 
is 7/10 the Bohr frequency. 



19 



FIG. 9. Free evolution of a two-dimensional gaussian wave packet. It has been obtained by solving the Schrodinger equation 
in the displayed rectangular region, i.e. horizontal axis x £ [—1,-1-1] and vertical axis y € [0,5], by imposing exact boundary 
conditions at finite distance. The countour plots of the probability density are shown for six equidistant times. The wave 
packet (t=l) was initially centered at x = y = 1 with an initial momentum v y /v x — 3/2. The figure shows the spreading of 
the wave packet in both directions and the crossing of the x = 1 boundary without any parasite reflections. 
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